## --------------------------------
# Main analysis                
## --------------------------------
# setwd("analysis/submission/")
source("help.r")

# load data
u_sample <- readRDS("u_sample.rds")

# list of covariates
ind_covs <- c("years_emp", "position", "educ", "college", "age", "income")
firm_covs <- c("industry", "industry_group", "manufacturing", "employee", "HQ_prefec", "TKY_OSK", "estab_year","foreign_owned", "foreign_sales", "product_sub", "indirect_import", "indirect_export", "capital",  "sales", "cptpp")
u_pretreat_covs = c("sanction", "sanction_impact","sanction_second")  

# parameters
treatment_repu_vec = c("u_treatment_firm_china_repu","u_treatment_firm_multiple_repu", "u_treatment_firm_us_repu")
treatment_vec = c("u_treatment_firm_china","u_treatment_firm_multiple", "u_treatment_firm_us")
control_vec = c("u_treatment_firm_multiple","u_treatment_firm_us", "u_treatment_control")


# convert to numeric outcomes
u_sample %>% 
  mutate_at(vars(starts_with("u_stop")), as.numeric) -> u_sample

# covariates
covariate_ind = c("years_emp","position","college","income")
covariate_firm = c("TKY_OSK","estab_year","foreign_owned","capital")
pretreat_covs = c("sanction", "sanction_impact", "sanction_second")

# china 
run_lm(
    data = u_sample %>% mutate(u_stop = as.numeric(u_stop)),
    treatment = "u_treatment_firm_china",
    control = "u_treatment_firm_multiple",
    outcome = "u_stop",
    covs = c(covariate_firm, covariate_ind,pretreat_covs),
    covs_include = TRUE) -> est1_china

# multiple
run_lm(
    data = u_sample %>% mutate(u_stop = as.numeric(u_stop)),
    treatment = "u_treatment_firm_multiple",
    control = "u_treatment_control",
    outcome = "u_stop",
    covs = c(covariate_firm, covariate_ind,pretreat_covs),
    covs_include = TRUE) -> est1_multi

# us
run_lm(
    data = u_sample %>% mutate(u_stop = as.numeric(u_stop)),
    treatment = "u_treatment_firm_us",
    control = "u_treatment_control",
    outcome = "u_stop",
    covs = c(covariate_firm, covariate_ind,pretreat_covs),
    covs_include = TRUE) -> est1_us

# plot numerical outcomes
bind_rows(
    broom::tidy(est1_us$fit) %>%
        mutate(names = "US withdrawal"),        
    broom::tidy(est1_multi$fit) %>%
        mutate(names = "Multiple withdrawal"),
    broom::tidy(est1_china$fit) %>%
        mutate(names = "China stays")) %>% 
    mutate(
        conf.low.95 = estimate - std.error * 1.96,
        conf.high.95 = estimate + std.error * 1.96,
        .width = 0.95,
        conf.low.90 = estimate - std.error * 1.64,
        conf.high.90 = estimate + std.error * 1.64,
        .width = 0.9) %>%
filter(startsWith(term, "u_treatment")) %>% 
mutate(
    cov = "With controls",
    names = factor(names, levels = c( "China stays", "US withdrawal", "Multiple withdrawal")),
    group = as_factor(c(rep("Peer Conformity", 2), "Market Competition"))) %>%
plot_ols_estimates(., "names") + 
my_theme() +
# change title text size to 8
theme(plot.title = element_text(size = 14)) +
ggtitle("Estimated Effects on Firm Managers'\n Support for Withdrawal from Russia") +
scale_x_continuous(breaks = c(-0.2, 0, 0.2)) +
theme(legend.position = "bottom") +
xlab("Changes in Support for Withdrawal from Russia") 
  


ggsave("plot/ukraine_H1-2_combined_ols.pdf", width = 6.2, height = 5.5)


# --------------------------------------------------------------
# with business in the market


# subgroup analysis

sub_us = u_sample %>% filter(trade_us_binary == 1) 

run_lm(
    data = sub_us,
    treatment = "u_treatment_firm_us",
    control = "u_treatment_control",
    outcome = "u_stop",
    covs = c(covariate_firm, covariate_ind,pretreat_covs),
    covs_include = TRUE) -> est1_us_sub

# wo business in US
sub_us2 = u_sample %>% filter(trade_us_binary == 0) 

run_lm(
    data = sub_us2,
    treatment = "u_treatment_firm_us",
    control = "u_treatment_control",
    outcome = "u_stop",
    covs = c(covariate_firm, covariate_ind,pretreat_covs),
    covs_include = TRUE) -> est1_us_sub2

# with business in other countries
sub_multi = u_sample %>% 
  filter(trade_multiple_binary == 1) 

run_lm(
    data = sub_multi,
    treatment = "u_treatment_firm_multiple",
    control = "u_treatment_control",
    outcome = "u_stop",
    covs = c(covariate_firm, covariate_ind,pretreat_covs),
    covs_include = TRUE) -> est1_multi_sub

# wo business in multiple countries
sub_multi2 = u_sample %>% 
  filter(trade_multiple_binary == 0) 

run_lm(
    data = sub_multi2,
    treatment = "u_treatment_firm_multiple",
    control = "u_treatment_control",
    outcome = "u_stop",
    covs = c(covariate_firm, covariate_ind,pretreat_covs),
    covs_include = TRUE) -> est1_multi_sub2

# with business in China
sub_china = u_sample %>% filter(trade_china_binary == 1) 

run_lm(
    data = sub_china,
    treatment = "u_treatment_firm_china",
    control = "u_treatment_firm_multiple",
    outcome = "u_stop",
    covs = c(covariate_firm, covariate_ind,pretreat_covs),
    covs_include = TRUE) -> est1_china_sub

# wo business in China
sub_china2 = u_sample %>% filter(trade_china_binary == 0) 

run_lm(
    data = sub_china2,
    treatment = "u_treatment_firm_china",
    control = "u_treatment_firm_multiple",
    outcome = "u_stop",
    covs = c(covariate_firm, covariate_ind,pretreat_covs),
    covs_include = TRUE) -> est1_china_sub2


# plot 
bind_rows(
    broom::tidy(est1_china_sub$fit) %>%
        mutate(legends = "China stays"),
    broom::tidy(est1_china_sub2$fit) %>%
        mutate(legends = "China stays"),        
    broom::tidy(est1_multi_sub$fit) %>%
        mutate(legends = "Multiple withdrawal"),
    broom::tidy(est1_multi_sub2$fit) %>%
        mutate(legends = "Multiple withdrawal"),
    broom::tidy(est1_us_sub$fit) %>%
        mutate(legends = "US withdrawal"),
    broom::tidy(est1_us_sub2$fit) %>%
        mutate(legends = "US withdrawal")
        ) %>%  
    mutate(
        conf.low.95 = estimate - std.error * 1.96,
        conf.high.95 = estimate + std.error * 1.96,
        .width = 0.95,
        conf.low.90 = estimate - std.error * 1.64,
        conf.high.90 = estimate + std.error * 1.64,
        .width = 0.9,
        legends = fct_relevel(legends, c(
            "US withdrawal", "Multiple withdrawal", "China stays"))) %>%
filter(startsWith(term, "u_treatment")) %>%
mutate(
    term = c(
        "china_business", "china_business_non",
        "multi_business", "multi_business_non",
        "us_business", "us_business_non"),
    cov = "With controls") -> data_fig3

# get p-value for interaction
u_sample %>% 
  filter(u_treatment_firm_us + u_treatment_control > 0) %>% 
  lm_robust(
    as.formula(
        paste("u_stop ~ u_treatment_firm_us*trade_us_binary +", 
        paste0(c(covariate_firm, covariate_ind, pretreat_covs), collapse = "+")) 
    ), 
    data = .,
    fixed_effects = ~ industry_group) %>%
    tidy() %>%
    filter(str_detect(term, ":")) %>%
    pull(p.value) -> p_value_us

u_sample %>%
    filter(u_treatment_firm_multiple + u_treatment_control > 0) %>%
    lm_robust(
        as.formula(
            paste("u_stop ~ u_treatment_firm_multiple*trade_multiple_binary +", 
            paste0(c(covariate_firm, covariate_ind, pretreat_covs), collapse = "+")) 
        ), 
        data = ., fixed_effects = ~ industry_group) %>%
        tidy() %>%
        filter(str_detect(term, ":")) %>%
        pull(p.value) -> p_value_multi

u_sample %>%
    filter(u_treatment_firm_china + u_treatment_firm_multiple > 0) %>%
    lm_robust(
        as.formula(
            paste("u_stop ~ u_treatment_firm_china*trade_china_binary +", 
            paste0(c(covariate_firm, covariate_ind, pretreat_covs), collapse = "+")) 
        ), 
        data = ., fixed_effects = ~ industry_group) %>%
        tidy() %>%
        filter(str_detect(term, ":")) %>%
        pull(p.value) -> p_value_china
    
# p-value table
pval_china <- tibble::tribble(
  ~group1, ~group2, ~p,     ~y.position,
  "china_business",    "china_business_non", round(p_value_china, 2), 0.1
)

pval_multi <- tibble::tribble(
  ~group1, ~group2, ~p,     ~y.position,
  "multi_business",    "multi_business_non", round(p_value_multi, 2), 0.1
)

pval_us <- tibble::tribble(
  ~group1, ~group2, ~p,     ~y.position,
  "us_business",    "us_business_non", round(p_value_us, 2), 0.1
)


# plot 
data_fig3 %>%
    filter(legends == "US withdrawal") %>% 
    ggplot(.) %>%
    plot_estimate_labels(
        ., pval_us, label = "p = {p} '", x_labs = "", y.position = 0.65) +
    geom_hline(
        yintercept =0, linetype = "dashed") +
    guides(color = FALSE) +
    scale_x_discrete(
        labels = c("w/o US\nbusiness", "US\nbusiness")) +
    scale_y_continuous(
        limits = c(-0.35, 0.7)) +
    ylab("Changes in Support for\nWithdrawal from Russia") +        
    ggtitle("US withdrawal") -> plot_us

data_fig3 %>% 
    filter(legends == "Multiple withdrawal") %>% 
    ggplot(.) %>%
    plot_estimate_labels(
        ., pval_multi, label = "p = {p} '", x_labs = "", y.position = 0.75) +
    geom_hline(
        yintercept =0, linetype = "dashed") +
    guides(color = FALSE) +
    scale_x_discrete(
        labels = c("w/o Foreign\nbusiness", "Foreign\nbusiness")) +
    ggtitle("Multiple withdrawal") +
    ylab("") +
    scale_y_continuous(
        limits = c(-0.2, 0.8)) -> plot_multi

data_fig3 %>% 
    filter(legends == "China stays")%>% 
    ggplot(.) %>%
    plot_estimate_labels(
        ., pval_china, label = "p = {p}", x_labs = "", y.position = 0.25) +
    geom_hline(
        yintercept = 0, linetype = "dashed") +
    guides(color = FALSE) +
    scale_x_discrete(
        labels = c("w/o Chinese\nbusiness", "Chinese\nbusiness")) +
    ggtitle("China stays") +
    ylab("") +
    scale_y_continuous(
        limits = c(-0.5, 0.3)) -> plot_china

plot_multi + plot_us  + plot_china -> plot_subgroup

plot_subgroup +
plot_annotation(
    title = 'Heterogeneous Effects on Support for Withdrawal from Russia by Market',
    theme = theme(plot.title = element_text(size = 20, hjust = 0.5, face="bold"))) 

ggsave("plot/ukraine_H1-2_subgroup_ols.pdf", width = 13, height = 5.5)


# ----------------------
# Reputation effect
# ----------------------

covariate_ind = c("years_emp","position","college","income")
covariate_firm = c("TKY_OSK","estab_year","foreign_owned","capital")
pretreat_covs = c("sanction", "sanction_impact", "sanction_second")

collect_models_repu = c()

for(i in 1:3){
  treatment = treatment_repu_vec[i]
  control = treatment_vec[i]

  run_lm(
    data = u_sample,
    treatment = treatment,
    control = control,
    outcome = "u_stop",
    covs = c(pretreat_covs, covariate_firm, covariate_ind),
    covs_include = TRUE) -> model_w
      
  collect_models_repu[[i]] = model_w
}

# get results
data_fig5 <- bind_rows(
  collect_models_repu[[3]] %>%
    get_data("Reputation effects\n(US withdraw branch)"),
  collect_models_repu[[2]] %>%
    get_data("Reputation effects\n(multiple withdraw branch)"),
  collect_models_repu[[1]] %>%
    get_data("Reputation effects\n(China stays branch)")
)

data_fig5 %>%
    mutate(
      label = factor(label, levels = c(
        "Reputation effects\n(China stays branch)",
        "Reputation effects\n(US withdraw branch)",
        "Reputation effects\n(multiple withdraw branch)")),
      group  = c(rep(
          "Peer Conformity", 2), "Market Competition") %>% as_factor,
      conf.low.95 = estimate - std.error * 1.96,
      conf.high.95 = estimate + std.error * 1.96,
      .width = 0.95,
      conf.low.90 = estimate - std.error * 1.64,
      conf.high.90 = estimate + std.error * 1.64,
      .width = 0.9) %>% 
plot_ols_estimates(., "label") +
my_theme() +
ggtitle("Reputation Effects on Support for\nWithdrawal from Russia") +
scale_x_continuous(breaks = c(-0.5, -0.4, -0.3, -0.2, -0.1, 0, 0.1, 0.2)) +
theme(legend.position = "bottom") +
xlab("Changes in Support for Withdrawal from Russia") 


ggsave("plot/ukraine_H4a-5a_ols.pdf",width = 7, height = 5.5)



# ----------------------
# reason for withdrawal
# ----------------------

u_sample %>% 
  filter(u_treatment_firm_china + u_treatment_firm_multiple + u_treatment_firm_us + u_treatment_control > 0) %>%
  dplyr::select(starts_with("mediate")) %>% colSums()/nrow(u_sample) -> med_sum

med_df = as_tibble(med_sum) %>%
  mutate(names = names(med_sum) %>% 
    str_replace(., "mediate_", "")) 

recode_names <- c(
  "repu_china" = "Reputation in Chinese market",
  "repu_others" = "Reputation among  other \nJapanese firms",
  "repu_west" = "Reputation in Western market",  
  "sanction" = "Secondary sanctions",
  "buss_risk" = "Business risks in Russia",
  "effective" = "Sanction effectiveness")

recode_category = c(
  "sanction" = "Non reputation",
  "repu_west" = "Reputation",
  "repu_china" = "Reputation",
  "repu_others" = "Reputation",
  "buss_risk" = "Non reputation",
  "effective" = "Non reputation")

med_df %>% 
  filter(!names %in% c("others", "guide_jp", "int_repu", "repu", "sanction_us", "sanction_china")) %>%
  mutate(
    category = recode_factor(names, !!!(recode_category)),
    names = recode_factor(names, !!!(recode_names))) %>%
  arrange((value)) %>%
  mutate(names = as_factor(names))  %>%
  ggplot(aes(x = names, y = value, fill = category)) + 
    coord_flip() +
    geom_bar(position="stack", stat="identity") +
    my_theme() +
    ylab("Proportion of Respondents") +
    xlab("") +
    scale_fill_manual(values = c("#273046","#FAD510"))


ggsave("plot/ukraine_mediator.pdf",width = 7, height = 5)



# ----------------------
# behavior outcomes
# ----------------------

# seek info -- public opinion
covariate_ind = c("years_emp","position","college","age")
covariate_firm = c("industry_group","employee","TKY_OSK","capital", "sales")
pretreat_covs = c("sanction", "sanction_impact", "sanction_second")

# run models
run_logit_all(
  "info_us_public", 
  "info_china_public", 
  "info_japan_public") -> out


# seek info -- gov policy
covariate_ind = c("years_emp","position","college","income")
covariate_firm = c("industry_group","employee","TKY_OSK","capital","sales")
pretreat_covs = c("sanction", "sanction_impact", "sanction_second")

# run models
run_logit_all(
  "info_us_sanction", 
  "info_china_statement", 
  "info_japan_support") -> out2

# US 
bind_rows(
  bind_rows(out$est1, out$est2) %>% 
    filter(term == "u_treatment_firm_us") %>%
    mutate(model = recode_factor(model, "US withdrawal" = "US public opinion")),
  bind_rows(out2$est1, out2$est2) %>% 
    filter(term == "u_treatment_firm_us") %>%
    mutate(model = recode_factor(model, "US withdrawal" = "US government policy"))
  ) %>% 
  mutate(
    cov = factor(model, levels = c( "US government policy", "US public opinion"))) %>% 
  plot_logit_estimates(., "model") +
  theme(legend.position= "none") +
  theme(axis.text.y = element_text(size = 15))+  
  ggtitle("US withdrawal") +
  xlab("") +
  scale_x_continuous(limits= c(-1, 1)) -> plot_info1

# China
bind_rows(
  bind_rows(out$est1, out$est2) %>% 
    filter(term == "u_treatment_firm_china") %>%
    mutate(model = recode_factor(model, "China stays" = "China's public opinion")),
  bind_rows(out2$est1, out2$est2) %>% 
    filter(term == "u_treatment_firm_china") %>%
    mutate(model = recode_factor(model, "China stays" = "Chinese government policy"))
  ) %>%
  mutate(
    cov = factor(model, levels = c("Chinese government policy", "China's public opinion"))) %>% 
  plot_logit_estimates(., "model") +
  ylab("") +
  theme(legend.position= "none") +
  ggtitle("China stays") +
  xlab("") +
  theme(axis.title.y=element_blank(),
        # axis.text.y=element_blank(),
        axis.ticks.y=element_blank(),
        plot.margin = margin(0.5, 0.5, 0.5, 0.5, "cm")) +
  scale_x_continuous(limits= c(-1, 1.3)) ->
  plot_info3

# combined plots
plot_info1 + plot_info3  -> plot_info

plot_info +
  plot_annotation(
    title = 'Heterogeneous Effects on Seeking Information',
    theme = theme(plot.title = element_text(size = 20, hjust = 0.7, face="bold"))) 

ggsave("plot/ukraine_H7-4_combined.pdf", width = 12, height = 5.5)

